	/****************************************************************************
	* 			Title: STYL Paper Figures
	*			Output:
	*			Date:
	*****************************************************************************/

	////1. DIRECTORY INFORMATION ////
	clear all
	set maxvar 30000
	use "data/clean/STYL_lr_panel_analysis.dta", clear

	////2. ADDITIONAL SET UP ////
*	set varabbrev off
*	set matsize 1000


	/////3.-  PREPARING APPENDIX DATA (BASELINE RESULTS AS ENDLINE WITH ROUND =  0)


	keep if round==1
	replace round=0


*Task1: Carry over baseline asb index as endline of round 0 to make the plot simpler.
	foreach var in uncens_fam_asb fam_asb cens_fam_asb {
		cap drop `var'_e
		ren  `var'_b `var'_e
	}


*Task2: Carry over and homologate pre-index components for the sake of re-standarization across all measures.
* We are using 6 (no domestic abuse Q at BL) of the rawest (no censored, std) versions of those questions
*asbhostil_b ~~ asbhostil_e
*disputes_all_b ~~ disputes_all_e
*carryweapon_b ~~ carryweapon_e
*drugssellever_b ~~ drugssellever_e
*stealnb_b      ~~ stealnb_e

	foreach var in asbhostil disputes_all carryweapon drugssellever stealnb  {
		cap drop `var'_e
		ren  `var'_b `var'_e
	}



keep 	partid           ///
		round            ///
        tpassonly        ///
        cashassonly      ///
        tpcashass        ///
        control          ///
        tp_strata_alt    ///
        cg_strata        ///
    uncens_fam_asb_e     ///
    fam_asb_e            ///
    cens_fam_asb_e       ///
    asbhostil_e            ///
    disputes_all_e         ///
    carryweapon_e          ///
    drugssellever_e        ///
    stealnb_e

	tempfile bl_data
	save  `bl_data'


	//// 4. ACTUALLY LOADING THE DATA ////

	use "data/clean/STYL_lr_panel_analysis.dta", clear


	//// 3. SURVEYSET, PROGRAMS AND GLOBALS ////

	    qui include dofiles/STYL_longrun_globals.do

	 	//c. SET SURVEY

	sort partid round
	svyset, strata(strata) psu(partid)


	    * Label Treatment arms
				qui la var tpassonly "Therapy Only"
				qui la var cashassonly "Cash Only"
				qui la var tpcashass "Therapy + Cash"

			/* Generate fixed effects dummies */
				qui tab tp_strata_alt, g(tp_dummy)
				qui tab cg_strata, g(cg_dummy)
				qui tab neighborhoodphase, g(nbhd_dummy)


	        cap drop __000001
	        cap drop __000000

*	egen patience_rstd_b = std(patient_game_real_b)
	append using `bl_data'




	foreach var in asbhostil disputes_all carryweapon drugssellever arrested stealnb {


*Top censoring across some measures across all rounds
			if inlist("`var'", "disputes_all", "stealnb") {

						sum `var'_e
						replace `var'_e = r(p99) if `var'_e > r(p99)
			}


*Normalize with respect to baseline mean.
		summ `var'_e if round==0
		gen `var'_wrt_blmean = (`var'_e-r(mean))/r(sd)


}





egen fam_asb_eb=rowmean(     ///
    asbhostil_wrt_blmean           ///
    disputes_all_wrt_blmean         ///
    carryweapon_wrt_blmean          ///
    drugssellever_wrt_blmean        ///
    stealnb_wrt_blmean              ///
 )


	* Second version where graphs match the tables (averages over rounds)
gen fam_asb_eb2=.
	replace fam_asb_eb2=fam_asb_e if round==0
	replace fam_asb_eb2=cens_fam_asb_stav if round==1 | round==2
	replace fam_asb_eb2=cens_fam_asb_ltav if round==5 | round==6
	replace fam_asb_eb2=cens_fam_asb_tyav if round==7 | round==8



****** Figure 1: ASB over time ******

************************************
****** Graph 1: ASB over time ******
*set scheme plotplainblind
*levels:

mat drop _all
* Treatment Arms - levels
mat define results_c = J(20,12, .)


reghdfe fam_asb_eb tpassonly cashassonly tpcashass  if round==0, abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb if round==0 & control==1
scalar ControlMean = r(mean)
dis ControlMean

mat results_c[1,1] = ControlMean
mat results_c[1,2] = r(lb)
mat results_c[1,3] = r(ub)

mat results_c[1,4] = _b[tpassonly]+ ControlMean
mat results_c[1,5] = _b[tpassonly]- (invttail(e(df_r),0.025)*_se[tpassonly]) + ControlMean
mat results_c[1,6] = _b[tpassonly]+ (invttail(e(df_r),0.025)*_se[tpassonly]) + ControlMean

mat results_c[1,7] = _b[cashassonly] + ControlMean
mat results_c[1,8] = _b[cashassonly]- (invttail(e(df_r),0.025)*_se[cashassonly]) + ControlMean
mat results_c[1,9] = _b[cashassonly]+ (invttail(e(df_r),0.025)*_se[cashassonly]) + ControlMean

mat results_c[1,10] = _b[tpcashass] + ControlMean
mat results_c[1,11] = _b[tpcashass]- (invttail(e(df_r),0.025)*_se[tpcashass]) + ControlMean
mat results_c[1,12] = _b[tpcashass]+ (invttail(e(df_r),0.025)*_se[tpcashass]) + ControlMean


reghdfe fam_asb_eb tpassonly cashassonly tpcashass $base if (round==1 | round==2), abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb   if (round==1| round==2) & control==1
scalar ControlMean = r(mean)
dis ControlMean

mat results_c[3,1] = ControlMean
mat results_c[3,2] = r(lb)
mat results_c[3,3] = r(ub)

mat results_c[3,4] = _b[tpassonly] + ControlMean
mat results_c[3,5] = _b[tpassonly]- (invttail(e(df_r),0.025)*_se[tpassonly]) + ControlMean
mat results_c[3,6] = _b[tpassonly]+ (invttail(e(df_r),0.025)*_se[tpassonly]) + ControlMean

mat results_c[3,7] = _b[cashassonly] + ControlMean
mat results_c[3,8] = _b[cashassonly]- (invttail(e(df_r),0.025)*_se[cashassonly]) + ControlMean
mat results_c[3,9] = _b[cashassonly]+ (invttail(e(df_r),0.025)*_se[cashassonly]) + ControlMean

mat results_c[3,10] = _b[tpcashass] + ControlMean
mat results_c[3,11] = _b[tpcashass]- (invttail(e(df_r),0.025)*_se[tpcashass]) + ControlMean
mat results_c[3,12] = _b[tpcashass]+ (invttail(e(df_r),0.025)*_se[tpcashass]) + ControlMean


reghdfe fam_asb_eb tpassonly cashassonly tpcashass $base if (round==5 | round==6), abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb   if (round==5| round==6) & control==1
scalar ControlMean = r(mean)
dis ControlMean

mat results_c[7,1] = ControlMean
mat results_c[7,2] = r(lb)
mat results_c[7,3] = r(ub)

mat results_c[7,4] = _b[tpassonly] + ControlMean
mat results_c[7,5] = _b[tpassonly]-(invttail(e(df_r),0.025)*_se[tpassonly]) + ControlMean
mat results_c[7,6] = _b[tpassonly]+(invttail(e(df_r),0.025)*_se[tpassonly]) + ControlMean

mat results_c[7,7] = _b[cashassonly] + ControlMean
mat results_c[7,8] = _b[cashassonly]-(invttail(e(df_r),0.025)*_se[cashassonly]) + ControlMean
mat results_c[7,9] = _b[cashassonly]+(invttail(e(df_r),0.025)*_se[cashassonly]) + ControlMean

mat results_c[7,10] = _b[tpcashass] + ControlMean
mat results_c[7,11] = _b[tpcashass]-(invttail(e(df_r),0.025)*_se[tpcashass]) + ControlMean
mat results_c[7,12] = _b[tpcashass]+(invttail(e(df_r),0.025)*_se[tpcashass]) + ControlMean


reghdfe fam_asb_eb tpassonly cashassonly tpcashass $base if (round==7 | round==8), abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb   if (round==7| round==8) & control==1
scalar ControlMean = r(mean)
dis ControlMean

mat results_c[20,1] = ControlMean
mat results_c[20,2] = r(lb)
mat results_c[20,3] = r(ub)

mat results_c[20,4] = _b[tpassonly] + ControlMean
mat results_c[20,5] = _b[tpassonly]-(invttail(e(df_r),0.025)*_se[tpassonly]) + ControlMean
mat results_c[20,6] = _b[tpassonly]+(invttail(e(df_r),0.025)*_se[tpassonly]) + ControlMean

mat results_c[20,7] = _b[cashassonly] + ControlMean
mat results_c[20,8] = _b[cashassonly]-(invttail(e(df_r),0.025)*_se[cashassonly]) + ControlMean
mat results_c[20,9] = _b[cashassonly]+(invttail(e(df_r),0.025)*_se[cashassonly]) + ControlMean

mat results_c[20,10] = _b[tpcashass] + ControlMean
mat results_c[20,11] = _b[tpcashass]-(invttail(e(df_r),0.025)*_se[tpcashass]) + ControlMean
mat results_c[20,12] = _b[tpcashass]+(invttail(e(df_r),0.025)*_se[tpcashass]) + ControlMean


matrix colname results_c = control lb ub coeftp lb ub coefcsh lb ub coeftpcsh lb ub


* ATE:
mat define results_d = J(20,12, .)


reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass if round==0, abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb2 if round==0 & control==1
scalar ControlMean = r(mean)
dis ControlMean

mat results_d[1,1] = ControlMean
mat results_d[1,2] = r(lb)
mat results_d[1,3] = r(ub)

mat results_d[1,4] = _b[tpassonly]
mat results_d[1,5] = _b[tpassonly]- (invttail(e(df_r),0.025)*_se[tpassonly])
mat results_d[1,6] = _b[tpassonly]+ (invttail(e(df_r),0.025)*_se[tpassonly])

mat results_d[1,7] = _b[cashassonly]
mat results_d[1,8] = _b[cashassonly]- (invttail(e(df_r),0.025)*_se[cashassonly])
mat results_d[1,9] = _b[cashassonly]+ (invttail(e(df_r),0.025)*_se[cashassonly])

mat results_d[1,10] = _b[tpcashass]
mat results_d[1,11] = _b[tpcashass]- (invttail(e(df_r),0.025)*_se[tpcashass])
mat results_d[1,12] = _b[tpcashass]+ (invttail(e(df_r),0.025)*_se[tpcashass])


reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==1 | round==2) , abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb2   if (round==1| round==2) & control==1
scalar ControlMean = r(mean)
dis ControlMean

mat results_d[3,1] = ControlMean
mat results_d[3,2] = r(lb)
mat results_d[3,3] = r(ub)

mat results_d[3,4] = _b[tpassonly]
mat results_d[3,5] = _b[tpassonly]- (invttail(e(df_r),0.025)*_se[tpassonly])
mat results_d[3,6] = _b[tpassonly]+ (invttail(e(df_r),0.025)*_se[tpassonly])

mat results_d[3,7] = _b[cashassonly]
mat results_d[3,8] = _b[cashassonly]- (invttail(e(df_r),0.025)*_se[cashassonly])
mat results_d[3,9] = _b[cashassonly]+ (invttail(e(df_r),0.025)*_se[cashassonly])

mat results_d[3,10] = _b[tpcashass]
mat results_d[3,11] = _b[tpcashass]- (invttail(e(df_r),0.025)*_se[tpcashass])
mat results_d[3,12] = _b[tpcashass]+ (invttail(e(df_r),0.025)*_se[tpcashass])


reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==5 | round==6), abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb2   if (round==5| round==6) & control==1
scalar ControlMean = r(mean)
dis ControlMean

mat results_d[7,1] = ControlMean
mat results_d[7,2] = r(lb)
mat results_d[7,3] = r(ub)

mat results_d[7,4] = _b[tpassonly]
mat results_d[7,5] = _b[tpassonly]-(invttail(e(df_r),0.025)*_se[tpassonly])
mat results_d[7,6] = _b[tpassonly]+(invttail(e(df_r),0.025)*_se[tpassonly])

mat results_d[7,7] = _b[cashassonly]
mat results_d[7,8] = _b[cashassonly]-(invttail(e(df_r),0.025)*_se[cashassonly])
mat results_d[7,9] = _b[cashassonly]+(invttail(e(df_r),0.025)*_se[cashassonly])

mat results_d[7,10] = _b[tpcashass]
mat results_d[7,11] = _b[tpcashass]-(invttail(e(df_r),0.025)*_se[tpcashass])
mat results_d[7,12] = _b[tpcashass]+(invttail(e(df_r),0.025)*_se[tpcashass])


reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==7 | round==8), abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb2   if (round==7| round==8) & control==1
scalar ControlMean = r(mean)
dis ControlMean

mat results_d[20,1] = ControlMean
mat results_d[20,2] = r(lb)
mat results_d[20,3] = r(ub)

mat results_d[20,4] = _b[tpassonly]
mat results_d[20,5] = _b[tpassonly]-(invttail(e(df_r),0.025)*_se[tpassonly])
mat results_d[20,6] = _b[tpassonly]+(invttail(e(df_r),0.025)*_se[tpassonly])

mat results_d[20,7] = _b[cashassonly]
mat results_d[20,8] = _b[cashassonly]-(invttail(e(df_r),0.025)*_se[cashassonly])
mat results_d[20,9] = _b[cashassonly]+(invttail(e(df_r),0.025)*_se[cashassonly])

mat results_d[20,10] = _b[tpcashass]
mat results_d[20,11] = _b[tpcashass]-(invttail(e(df_r),0.025)*_se[tpcashass])
mat results_d[20,12] = _b[tpcashass]+(invttail(e(df_r),0.025)*_se[tpcashass])


matrix colname results_d = control lb ub coeftp lb ub coefcsh lb ub coeftpcsh lb ub
*matrix rowname results_b = fam_asb_b fam_asb_st fam_asb_lt fam_asb_ty


graph set window fontface "Times New Roman"


coefplot (matrix(results_c[.,1]), noci  label("Control") recast(connected) msymbol(X) msize(medium) mcolor(red) lcolor(red%50)  lwidth(medthick)) ///
     (matrix(results_c[.,4]), noci label("Therapy") recast(connected) msymbol(X) msize(medium)  mcolor(lavender) lcolor(lavender%50)) ///
     (matrix(results_c[.,7]), noci label("Cash") recast(connected) msymbol(X) msize(medium)  mcolor(navy) lcolor(navy%50)) ///
     (matrix(results_c[.,10]), noci label("Therapy and Cash") recast(connected) msymbol(X) msize(medium)  mcolor(green) lcolor(green%50) ), ///
      vertical graphregion(fcolor(white)) legend(row(1) size(small)) xlabel(1 "Baseline" 3 "1 Month" 7 "1 Year" 20 "10 years", labs(small) angle(45)) ///
      ylabel(-0.4(0.2)0.2) yscale(range(-0.3, 0.1) noextend)  xsize(5) ysize(5) title("a) Levels, All Treatment Arms", size(medsmall)) legend(off)  ///
      text(-0.25 17 "Control",  color(red) place(n) size(small)) text(-0.2 10 "Cash",  color(navy) place(n) size(small)) ///
      text(-0.3 10 "Therapy", color(lavender) place(n) size(small)) text(-0.4 14 "Therapy and Cash", color(green) place(n)size(small))  ytitle("Standard deviations")
graph save graph1d.gph, replace

coefplot  (matrix(results_d[.,4]), ci((results_d[.,5] results_d[.,6])) label("Therapy") recast(connected) msymbol(X) msize(medium) mcolor(lavender) lcolor(lavender%50) ///
          ciopts(recast(rarea) color(lavender%25))), ///
      vertical yline(0,lcolor(red)) graphregion(fcolor(white)) legend(row(1) size(small)) xlabel(1 "Baseline" 3 "1 Month" 7 "1 Year" 20 "10 years", labs(small) angle(45)) ///
      ylabel(-0.4(0.2)0.2) yscale(range(-0.3, 0.1)noextend)  xsize(5) ysize(5) title("b) ATE of Therapy Only", size(medsmall)) legend(off)
graph save graph2d.gph, replace

coefplot  (matrix(results_d[.,7]), ci((results_d[.,8] results_d[.,9])) label("Cash") recast(connected) msymbol(X) msize(medium) mcolor(navy) lcolor(navy%50) ///
          ciopts(recast(rarea) color(navy%25))), ///
      vertical yline(0,lcolor(red)) graphregion(fcolor(white)) legend(row(1) size(small)) xlabel(1 "Baseline" 3 "1 Month" 7 "1 Year" 20 "10 years", labs(small) angle(45)) ///
      ylabel(-0.4(0.2)0.2) yscale(range(-0.3, 0.1)noextend)  xsize(5) ysize(5) title("c) ATE of Cash Only", size(medsmall)) legend(off)  ytitle("Standard deviations")
graph save graph3d.gph, replace

coefplot  (matrix(results_d[.,10]), ci((results_d[.,11] results_d[.,12])) label("Therapy + Cash") recast(connected) msymbol(X) msize(medium)  mcolor(green) lcolor(green%50) ///
          ciopts(recast(rarea) color(green%25))), ///
      vertical yline(0,lcolor(red)) graphregion(fcolor(white)) legend(row(1) size(small)) xlabel(1 "Baseline" 3 "1 Month" 7 "1 Year" 20 "10 years", labs(small) angle(45)) ///
      ylabel(-0.4(0.2)0.2) yscale(range(-0.3, 0.1)noextend)  xsize(5) ysize(5) title("d) ATE of Therapy + Cash", size(medsmall)) legend(off)
graph save graph4d.gph, replace

*grc1leg graph1b.gph graph2b.gph graph3b.gph graph4b.gph, col(2) xsize(5) ysize(4.5) ycommon xcommon
graph combine graph1d.gph graph2d.gph graph3d.gph graph4d.gph, col(2) xsize(5) ysize(4.5)  graphregion(color(white) lwidth(thick))
graph export "outfiles/figures//Figure1_Level_ATE.pdf", replace


erase graph1d.gph
erase graph2d.gph
erase graph3d.gph
erase graph4d.gph



* *******************************
****** Figure 2: ASB over time by baseline asb
*Het analysis Bruno: Check bl_asb_b differences.


astile quart_asb_b = fam_asb_b if round == 1, nquantiles(4)
bysort partid: egen asb_b_Q =max(quart_asb_b)


gen asb_b_Q4 = (asb_b_Q == 4)


* HIGH ASB
mat define results_d = J(20,12, .)

* reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass if round==0, abs(tp_strata_alt cg_strata)
reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass if round==0 & asb_b_Q4==1, abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb2 if round==0 & control==1 & asb_b_Q4==1
scalar ControlMean = r(mean)
dis ControlMean

mat results_d[1,1] = ControlMean
mat results_d[1,2] = r(lb)
mat results_d[1,3] = r(ub)

mat results_d[1,4] = _b[tpassonly]
mat results_d[1,5] = _b[tpassonly]- (invttail(e(df_r),0.025)*_se[tpassonly])
mat results_d[1,6] = _b[tpassonly]+ (invttail(e(df_r),0.025)*_se[tpassonly])

mat results_d[1,7] = _b[cashassonly]
mat results_d[1,8] = _b[cashassonly]- (invttail(e(df_r),0.025)*_se[cashassonly])
mat results_d[1,9] = _b[cashassonly]+ (invttail(e(df_r),0.025)*_se[cashassonly])

mat results_d[1,10] = _b[tpcashass]
mat results_d[1,11] = _b[tpcashass]- (invttail(e(df_r),0.025)*_se[tpcashass])
mat results_d[1,12] = _b[tpcashass]+ (invttail(e(df_r),0.025)*_se[tpcashass])


* reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==1 | round==2) , abs(tp_strata_alt cg_strata)
reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==1 | round==2) & asb_b_Q4==1, abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb2   if (round==1| round==2) & control==1 & asb_b_Q4==1
scalar ControlMean = r(mean)
dis ControlMean

mat results_d[3,1] = ControlMean
mat results_d[3,2] = r(lb)
mat results_d[3,3] = r(ub)

mat results_d[3,4] = _b[tpassonly]
mat results_d[3,5] = _b[tpassonly]- (invttail(e(df_r),0.025)*_se[tpassonly])
mat results_d[3,6] = _b[tpassonly]+ (invttail(e(df_r),0.025)*_se[tpassonly])

mat results_d[3,7] = _b[cashassonly]
mat results_d[3,8] = _b[cashassonly]- (invttail(e(df_r),0.025)*_se[cashassonly])
mat results_d[3,9] = _b[cashassonly]+ (invttail(e(df_r),0.025)*_se[cashassonly])

mat results_d[3,10] = _b[tpcashass]
mat results_d[3,11] = _b[tpcashass]- (invttail(e(df_r),0.025)*_se[tpcashass])
mat results_d[3,12] = _b[tpcashass]+ (invttail(e(df_r),0.025)*_se[tpcashass])


*reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==5 | round==6), abs(tp_strata_alt cg_strata)
reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==5 | round==6) & asb_b_Q4==1, abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb2   if (round==5| round==6) & control==1 & asb_b_Q4==1
scalar ControlMean = r(mean)
dis ControlMean

mat results_d[7,1] = ControlMean
mat results_d[7,2] = r(lb)
mat results_d[7,3] = r(ub)

mat results_d[7,4] = _b[tpassonly]
mat results_d[7,5] = _b[tpassonly]-(invttail(e(df_r),0.025)*_se[tpassonly])
mat results_d[7,6] = _b[tpassonly]+(invttail(e(df_r),0.025)*_se[tpassonly])

mat results_d[7,7] = _b[cashassonly]
mat results_d[7,8] = _b[cashassonly]-(invttail(e(df_r),0.025)*_se[cashassonly])
mat results_d[7,9] = _b[cashassonly]+(invttail(e(df_r),0.025)*_se[cashassonly])

mat results_d[7,10] = _b[tpcashass]
mat results_d[7,11] = _b[tpcashass]-(invttail(e(df_r),0.025)*_se[tpcashass])
mat results_d[7,12] = _b[tpcashass]+(invttail(e(df_r),0.025)*_se[tpcashass])


* reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==7 | round==8), abs(tp_strata_alt cg_strata)
reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==7 | round==8) & asb_b_Q4==1, abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb2   if (round==7| round==8) & control==1 & asb_b_Q4==1
scalar ControlMean = r(mean)
dis ControlMean

mat results_d[20,1] = ControlMean
mat results_d[20,2] = r(lb)
mat results_d[20,3] = r(ub)

mat results_d[20,4] = _b[tpassonly]
mat results_d[20,5] = _b[tpassonly]-(invttail(e(df_r),0.025)*_se[tpassonly])
mat results_d[20,6] = _b[tpassonly]+(invttail(e(df_r),0.025)*_se[tpassonly])

mat results_d[20,7] = _b[cashassonly]
mat results_d[20,8] = _b[cashassonly]-(invttail(e(df_r),0.025)*_se[cashassonly])
mat results_d[20,9] = _b[cashassonly]+(invttail(e(df_r),0.025)*_se[cashassonly])

mat results_d[20,10] = _b[tpcashass]
mat results_d[20,11] = _b[tpcashass]-(invttail(e(df_r),0.025)*_se[tpcashass])
mat results_d[20,12] = _b[tpcashass]+(invttail(e(df_r),0.025)*_se[tpcashass])


matrix colname results_d = control lb ub coeftp lb ub coefcsh lb ub coeftpcsh lb ub
*matrix rowname results_b = fam_asb_b fam_asb_st fam_asb_lt fam_asb_ty



* LOW ASB
mat define results_c = J(20,12, .)

* reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass if round==0, abs(tp_strata_alt cg_strata)
reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass if round==0 & asb_b_Q4==0, abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb2 if round==0 & control==1 & asb_b_Q4==0
scalar ControlMean = r(mean)
dis ControlMean

mat results_c[1,1] = ControlMean
mat results_c[1,2] = r(lb)
mat results_c[1,3] = r(ub)

mat results_c[1,4] = _b[tpassonly]
mat results_c[1,5] = _b[tpassonly]- (invttail(e(df_r),0.025)*_se[tpassonly])
mat results_c[1,6] = _b[tpassonly]+ (invttail(e(df_r),0.025)*_se[tpassonly])

mat results_c[1,7] = _b[cashassonly]
mat results_c[1,8] = _b[cashassonly]- (invttail(e(df_r),0.025)*_se[cashassonly])
mat results_c[1,9] = _b[cashassonly]+ (invttail(e(df_r),0.025)*_se[cashassonly])

mat results_c[1,10] = _b[tpcashass]
mat results_c[1,11] = _b[tpcashass]- (invttail(e(df_r),0.025)*_se[tpcashass])
mat results_c[1,12] = _b[tpcashass]+ (invttail(e(df_r),0.025)*_se[tpcashass])


* reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==1 | round==2) , abs(tp_strata_alt cg_strata)
reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==1 | round==2) & asb_b_Q4==0, abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb2   if (round==1| round==2) & control==1 & asb_b_Q4==0
scalar ControlMean = r(mean)
dis ControlMean

mat results_c[3,1] = ControlMean
mat results_c[3,2] = r(lb)
mat results_c[3,3] = r(ub)

mat results_c[3,4] = _b[tpassonly]
mat results_c[3,5] = _b[tpassonly]- (invttail(e(df_r),0.025)*_se[tpassonly])
mat results_c[3,6] = _b[tpassonly]+ (invttail(e(df_r),0.025)*_se[tpassonly])

mat results_c[3,7] = _b[cashassonly]
mat results_c[3,8] = _b[cashassonly]- (invttail(e(df_r),0.025)*_se[cashassonly])
mat results_c[3,9] = _b[cashassonly]+ (invttail(e(df_r),0.025)*_se[cashassonly])

mat results_c[3,10] = _b[tpcashass]
mat results_c[3,11] = _b[tpcashass]- (invttail(e(df_r),0.025)*_se[tpcashass])
mat results_c[3,12] = _b[tpcashass]+ (invttail(e(df_r),0.025)*_se[tpcashass])


*reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==5 | round==6), abs(tp_strata_alt cg_strata)
reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==5 | round==6) & asb_b_Q4==0, abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb2   if (round==5| round==6) & control==1 & asb_b_Q4==0
scalar ControlMean = r(mean)
dis ControlMean

mat results_c[7,1] = ControlMean
mat results_c[7,2] = r(lb)
mat results_c[7,3] = r(ub)

mat results_c[7,4] = _b[tpassonly]
mat results_c[7,5] = _b[tpassonly]-(invttail(e(df_r),0.025)*_se[tpassonly])
mat results_c[7,6] = _b[tpassonly]+(invttail(e(df_r),0.025)*_se[tpassonly])

mat results_c[7,7] = _b[cashassonly]
mat results_c[7,8] = _b[cashassonly]-(invttail(e(df_r),0.025)*_se[cashassonly])
mat results_c[7,9] = _b[cashassonly]+(invttail(e(df_r),0.025)*_se[cashassonly])

mat results_c[7,10] = _b[tpcashass]
mat results_c[7,11] = _b[tpcashass]-(invttail(e(df_r),0.025)*_se[tpcashass])
mat results_c[7,12] = _b[tpcashass]+(invttail(e(df_r),0.025)*_se[tpcashass])


* reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==7 | round==8), abs(tp_strata_alt cg_strata)
reghdfe fam_asb_eb2 tpassonly cashassonly tpcashass $base if (round==7 | round==8) & asb_b_Q4==0, abs(tp_strata_alt cg_strata)

ci mean fam_asb_eb2   if (round==7| round==8) & control==1 & asb_b_Q4==0
scalar ControlMean = r(mean)
dis ControlMean

mat results_c[20,1] = ControlMean
mat results_c[20,2] = r(lb)
mat results_c[20,3] = r(ub)

mat results_c[20,4] = _b[tpassonly]
mat results_c[20,5] = _b[tpassonly]-(invttail(e(df_r),0.025)*_se[tpassonly])
mat results_c[20,6] = _b[tpassonly]+(invttail(e(df_r),0.025)*_se[tpassonly])

mat results_c[20,7] = _b[cashassonly]
mat results_c[20,8] = _b[cashassonly]-(invttail(e(df_r),0.025)*_se[cashassonly])
mat results_c[20,9] = _b[cashassonly]+(invttail(e(df_r),0.025)*_se[cashassonly])

mat results_c[20,10] = _b[tpcashass]
mat results_c[20,11] = _b[tpcashass]-(invttail(e(df_r),0.025)*_se[tpcashass])
mat results_c[20,12] = _b[tpcashass]+(invttail(e(df_r),0.025)*_se[tpcashass])


matrix colname results_c = control lb ub coeftp lb ub coefcsh lb ub coeftpcsh lb ub


* Graph
coefplot (matrix(results_c[.,4]), noci label("Low ASB") recast(connected) msymbol(triangle) msize(medium) mcolor(gray) lcolor(gray%50)) ///
         (matrix(results_d[.,4]), ci((results_d[.,5] results_d[.,6]))  label("High ASB") recast(connected) msymbol(X)  msize(large) mcolor(lavender) lcolor(lavender%50) ///
          ciopts(recast(rarea) color(lavender%25))), vertical yline(0,lcolor(red)) graphregion(fcolor(white)) legend(row(1) size(small)) ///
           xlabel(1 "Baseline" 3 "1 Month" 7 "1 Year" 20 "10 years", labs(medium) angle(45)) ylabel(-1.0(0.5)1.1) yscale(range(-1.0, 1.1)noextend) ///
            xsize(5) ysize(5) title("a) ATE of Therapy Only", size(large))  ytitle("Standard deviations")
graph save graph2a.gph, replace

coefplot (matrix(results_c[.,7]), noci label("Low ASB") recast(connected) msymbol(triangle) mcolor(gray) lcolor(gray%50)) ///
         (matrix(results_d[.,7]), ci((results_d[.,8] results_d[.,9])) label("High ASB") recast(connected) msymbol(X) msize(large)  mcolor(navy) lcolor(navy%50) ///
          ciopts(recast(rarea) color(navy%25))),  vertical yline(0,lcolor(red)) graphregion(fcolor(white)) legend(row(1) size(small)) ///
           xlabel(1 "Baseline" 3 "1 Month" 7 "1 Year" 20 "10 years", labs(medium) angle(45)) ylabel(-1.0(0.5)1.1) yscale(range(-1.0, 1.1)noextend) ///
            xsize(5) ysize(5) title("b) ATE of Cash Only", size(large)) ytitle("Standard deviations")
graph save graph2b.gph, replace

coefplot (matrix(results_c[.,10]), noci label("Low ASB") recast(connected) msymbol(triangle) mcolor(gray) lcolor(gray%50)) ///
         (matrix(results_d[.,10]), ci((results_d[.,11] results_d[.,12])) label("High ASB") recast(connected) msymbol(X) msize(large)  mcolor(green) lcolor(green%50) ///
          ciopts(recast(rarea) color(green%25))),  vertical yline(0,lcolor(red)) graphregion(fcolor(white)) legend(row(1) size(small)) ///
           xlabel(1 "Baseline" 3 "1 Month" 7 "1 Year" 20 "10 years", labs(medium) angle(45)) ylabel(-1.0(0.5)1.1) yscale(range(-1.0, 1.1)noextend) ///
            xsize(5) ysize(5) title("c) ATE of Therapy + Cash", size(large)) ytitle("Standard deviations")
graph save graph2c.gph, replace

*grc1leg2       graph2a.gph graph2b.gph graph2c.gph, col(1) xsize(3) ysize(6)
graph combine graph2a.gph graph2b.gph graph2c.gph , col(1) xsize(3) ysize(6)  graphregion(color(white) lwidth(thick))

graph export "outfiles/figures//Figure2_ATE_het_ASB.pdf", replace
erase graph2a.gph
erase graph2b.gph
erase graph2c.gph
